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THE EFFECT OF GAS COOLING ON THE SHAPES OF DARK MATTER HALOS 
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ABSTRACT 

We analyze the effect of dissipation on the shapes of dark matter (DM) halos using high-resolution cosmological 
gasdynamics simulations of clusters and galaxies in the ACDM cosmology. We find that halos formed in simula- 
tions with gas cooling are significantly more spherical than corresponding halos formed in adiabatic simulations. 
Gas cooling results in an average increase of the principle axis ratios of halos by ~ 0.2 - 0.4 in the inner regions. 
The systematic difference decreases slowly with radius but persists almost to the virial radius. We argue that the 
differences in simulations with and without cooling arise both during periods of quiescent evolution, when gas 
cools and condenses toward the center, and during major mergers. We perform a series of high-resolution A-body 
simulations to study the shapes of remnants in major mergers of DM halos and halos with embedded stellar disks. 
In the DM halo-only mergers, the shape of the remnants depends only on the orbital angular momentum of the 
encounter and not on the internal structure of the halos. However, significant shape changes in the DM distribution 
may result if stellar disks are included. In this case the shape of the DM halos is correlated with the morphology 
of the stellar remnants. 

Subject headings: cosmology: theory — dark matter — galaxies: halos — halos: shapes — halos: structure — 
methods: numerical 



1. INTRODUCTION 

Triaxial dark matter (DM) halos are a generic prediction of 
the hierarchical, cold dark matter (CDM) model of structure 
formation. Therefore, observational probes of halo shapes are 
fundamental tests of this model. Recently, such comparisons 
have received attention as the accuracy and number of observa- 
tional probes of halo shapes increase for the Milky Way (MW) 
(see Merrifield 2003, for a recent review), other galaxies (e.g., 
Buote & Canizares 1998; Buote et al. 2002; Hoekstra et al. 
2004), and galaxy clusters (e.g., Kolokotronis et al. 2001). 

Dubinski & Carlberg (1991) first pointed out that observed 
elliptical galaxies are systematically more spherical than simu- 
lated CDM halos. It was recently argued that the coherence of 
the tidal stream of the Sagittarius dwarf spheroidal galaxy in- 
dicates that the inner halo of the MW is nearly spherical, with 
a minor- to-major axis ratio c/a > 0.8 (Ibata et al. 2001; Ma- 
jewski et al. 2003, but see Mayer et al. 2002, Helmi 2004, 
and Martmez-Delgado et al. 2004). This differs considerably 
from an average predicted ratio of (c/a) ~ 0.6-0.7 for MW- 
sized halos formed in dissipationless cosmological simulations 
(Barnes & Efstathiou 1987; Frenk et al. 1988; Dubinski & Carl- 
berg 1991; Warren et al. 1992; Cole & Lacey 1996; Bullock 
2002; Jing & Suto 2002). Although the axis ratios of simulated 
halos exhibit a broad distribution, with a Gaussian rms disper- 
sion of ~ 0.1, we would have to accept the MW halo as an 
outlier. 

Utilizing halo shapes as a test of the CDM paradigm requires 
improving both observations and theoretical predictions. Most 
of the numerical work published to date does not account for 
the effects of gas cooling on halo shapes. Notable exceptions 
are the studies by Katz & Gunn (1991), Evrard et al. (1994), 
and Tissera & Dommguez-Tenreiro (1998), who find that halos 
in dissipational simulations are systematically more spherical 
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than corresponding halos in dissipationless runs. However, the 
halos in these studies were resolved with only a few hundred 
to a few thousand particles and the shapes of the inner regions 
of halos could not be studied reliably. In addition, these stud- 
ies did not address the processes by which cooling affects halo 
shapes. Recently, Floor et al. (2003) used gasdynamics sim- 
ulations of galaxy clusters with and without radiative cooling 
to investigate the evolution of eccentricity. Their study showed 
very slow eccentricity evolution in simulated clusters. How- 
ever, they focused on the eccentricity at large radii where, as 
we demonstrate, the effect of cooling is small. 

Dubinski (1994) studied the effects of dissipation on halo 
shapes numerically by adiabatically growing central mass con- 
centrations in an initially triaxial A-body halo. This study showed 
that the growth of a central condensation results in a consid- 
erably rounder halo shape. These experiments, although sug- 
gestive, were not generally representative of the cosmological 
framework of hierarchical halo assembly via multiple mergers. 

In this Letter we revisit the effects of dissipation on halo 
shapes, combining self-consistent cosmological simulations with 
a resolution up to 2 orders of magnitude higher than in previous 
studies, with high-resolution controlled A-body simulations of 
binary mergers of disk galaxies and DM halos whose structure 
is based on the currently favored galaxy formation paradigm. 

2. NUMERICAL SIMULATIONS 

We analyze high-resolution cosmological simulations of eight 
group- and cluster-size systems and one galaxy-size system in a 
flat ACDM model: Vl m = 1 - A = 0.3, il h = 0.043, h = 0.7 and 
og = 0.9. The simulations were performed with the Adaptive 
Refinement Tree A-body+gasdynamics code (Kravtsov 1999; 
Kravtsov et al. 2002), an Eulerian code that uses adaptive re- 
finement in space and time and (nonadaptive) refinement in 
mass to achieve the high dynamic range required to resolve 
cores of halos formed in cosmological simulations. 

The cluster simulations have a peak resolution of « 2.44 /i" 1 
kpc and DM particle mass of 2.7 x 10 8 /i _I M Q with only a re- 
gion of ~ 10/i -1 Mpc around each cluster adaptively refined. 
We analyze each cluster at an epoch near z = 0, when it appears 
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most relaxed. This minimizes the noise introduced by substruc- 
ture, which may cause transient changes in the estimates of axis 
ratios at particular radii. For one cluster, we also analyze its 
progenitors at z = 1 and z = 2. The virial masses 4 of the clus- 
ters range from i=s 10 13 to 3 x 10 14 /i _1 M Q . The galaxy forma- 
tion simulation follows the early (z > 4) evolution of a galaxy 
that becomes an MW-size object at z = in a periodic box of 
6 h~ l Mpc. At z = 4, the galaxy already contains a large fraction 
of its final mass: w2x 10" h M©. The DM particle mass is 
9.18 x 10 5 /i _1 Mq and the peak resolution of the simulation is 
183 /i~ 1 comoving pc. This simulation is presented in Kravtsov 
(2003), where more details can be found. 

For each object, we analyze two sets of simulations started 
from the same initial conditions but including different physi- 
cal processes. The first set of simulations follows the dynam- 
ics of gas "adiabatically," i.e., without radiative cooling. For 
the cluster halos we also analyze purely dissipationless simula- 
tions without gas the results of which agree well with the adia- 
batic gasdynamics runs. The second set of simulations includes 
star formation, metal enrichment and thermal supernovae feed- 
back, metallicity- and density-dependent cooling, and UV heat- 
ing due to a cosmological ionizing background. We use these 
cosmological simulations to study the effect of dissipation on 
the shapes of halos in the next section. 

To explore the mechanisms responsible for the effect of cool- 
ing on halo shapes, we performed controlled iV-body merger 
simulations of equal-mass pure DM halos and multicomponent 
galaxies varying only the internal mass distributions (Kazantzidis 
et al., in preparation). The DM halos followed the Navarro et al. 
(1996, hereafter NFW) density profile and were initially spher- 
ical with isotropic velocity dispersion tensors (Kazantzidis et al. 
2004). First, we simulated several mergers of systems on parabolic 
orbits with pericenters that were 20% of the halo virial radii, 
typical of cosmological mergers (Khochfar & Burkert 2003). 
In one set of experiments, we merged DM halos with embed- 
ded stellar disks (Hernquist 1993; Mo et al. 1998). In these 
experiments we used randomly inclined and coplanar disk con- 
figurations and adopted halo parameters from the MW model 
Al of Klypin et al. (2002). Specifically, each halo had a virial 
mass of M v j r = 7 x IO'^^Mq, a concentration parameter of 
c = 12, and a dimensionless spin parameter of A = 0.031. The 
halos were adiabatically contracted to respond to the growth of 
the stellar disk. The mass and thickness of the stellar disk were 
M d = 0.04M vir and zo = 0. 17? d , respectively, and 7? d = 2.7 h~ x kpc 
was the disk scale length. In another experiment, we studied 
the merging of identical NFW halos with no disk component. 
In a third experiment, we studied the merging of "contracted 
halos" with steepened inner density profiles (p ~ r ). These 
profiles were set to match the combined initial spherical density 
profile of the adiabatically contracted halos and stellar disks in 
the halo+disk mergers. Last, we repeated the halo+disk and 
pure NFW halo mergers for radial orbits to study the effect of 
the impact parameter on the structure of the remnant (Moore 
et al. 2003). The initial separation of the halo centers was twice 
their virial radii. Their initial relative velocity was determined 
from the corresponding Keplerian orbit of two point masses for 
the parabolic mergers, while for the radial mergers it was set to 
the velocity of the circular orbit about their common center of 
mass. 

The merger simulations were performed with the PKDGRAV 

4 We define the virial radius as the radius enclosing an average density of 180 
times the mean density of the universe at the epoch of analysis. 
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FIG. 1. — Top panels: Axis ratios cja vs. b/a at 10% and 30% of the virial 
radius [using bins of A(r/r v ; r ) = 0.1] for halos in cosmological simulations. 
Open symbols show adiabatic simulations and filled symbols correspond to 
simulations with gas cooling and star formation. The galaxy simulation is 
shown by a square and the clusters by triangles. Bottom panel: Average dif- 
ference between axis ratios in cooling and adiabatic simulations as a function 
of radius. The error bars show the \o scatter about the mean value of A(b/a) 
in each bin. The scatter in A(c/a) is similar. 

(Stadel 2001), multistepping, parallel A^-body tree code. We 
used N = 2 x 10 5 DM particles and N = 2x 10 4 collisionless 
stellar particles and employed a gravitational softening of 0.7 
and 0.3 h~ l kpc, respectively. Our results did not change when 
we used a factor of 10 more particles and half the softening 
indicating the achievement of numerical convergence. 

3. RESULTS 

For each object, we calculate principle axis ratios s = b/a and 
q = c/a (a > b > c), from the eigenvalues of a modified iner- 
tia tensor (e.g., Dubinski & Carlberg 1991): Ijj = J2 a x ? x f Aa> 
where xf is the i coordinate of the ath particle, r 2 a = (yf) 2 + 
(yf I s ) 2 + (yf I q) 2 > an d yf are coordinates with respect to the 
principle axes. We use an iterative algorithm starting with a 
spherical configuration (a = b = c) and use the results of the pre- 
vious iteration to define the principle axes of the next iteration 
until the results converge to a fractional difference of 10~ 3 . We 
compute axis ratios "differentially" using only particles within 
finite bins of r, and label each bin by the mean value of r of the 
particles in the bin. We use differential measurements because 
the axis ratios calculated at different values of r are almost in- 
dependent. We find that in the standard practice of measuring 
shapes cumulatively and weighting each particle contribution 
to Ijj by r~ 2 , axis ratios at large r are quite sensitive to the dis- 
tribution of particles in the inner regions of the halo. The dif- 
ferences between the differential and cumulative estimates of 
axis ratios can be as large as a few tenths at large radii. How- 
ever, we present shape measurements in differential radial bins. 
Within each bin, different definitions of /y result in very small 
differences. 

Figure 1 shows the principle axis ratios of the dark matter at 
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FIG. 2. — Density maps for the galaxy-size (top panels) and a cluster-size 
(bottom panels) halo projected onto the plane of intermediate and major axis. 
The particles are color-coded on a logarithmic scale with brighter colors in 
higher density regions. The local density is calculated using an SPH smoothing 
kernel of 64 particles to suppress small substructure and to better show the 
shape of the diffuse DM in the main halo. Halos in the cooling runs (right 
panels) are significantly rounder than their adiabatic counterparts (left panels). 



10% and 30% of the virial radius for all cosmological simula- 
tions. The purely dissipationless simulations without gas pro- 
duce results similar to the adiabatic runs and are omitted for 
clarity. Note that the halos in adiabatic simulations tend to be 
slightly more elongated closer to their centers. Interestingly, 
there is a correlation of c/a and b/a. The effect of cooling 
changes axis ratios roughly along the correlation in the c fa- 
b/a plane by up to ~ 0.4 at 0.1r v i r and by ~ 0.2-0.3 at 0.3r v ; r . 
The significance of the change in shape can be seen in Figure 2, 
where we show density maps of the DM distribution in the sim- 
ulated galaxy (top panels) and a cluster (bottom panels) at z = 4 
and z = 0.4, respectively. The difference between the adiabatic 
and cooling simulations is visually striking, with the DM halos 
significantly rounder in the cooling runs. 

The radial dependence of the effect can be seen in the bot- 
tom panel of Figure 1, which shows the average difference be- 
tween axis ratios in the cooling and adiabatic runs [A(b/a) = 
(b/a) coo \-(b/a) a dvdb] as a function of radius. Although the ef- 
fect decreases with radius, the axis ratios in simulations with 
cooling are systematically larger than in the adiabatic and dissi- 
pationless cases, even at r « 0.8r v ; r . Figure 3 shows individual 
radial profiles of the minor-to-major axis ratio c/a in the galaxy 
and cluster simulations, as well as in representative controlled 
merger experiments. It is interesting to compare cluster and 
galaxy halos because their baryon distributions have different 
morphologies. In the galaxy, most (~ 90%) of the baryons lie in 
a highly flattened gaseous disk, while in the cluster simulations 
most of the baryons are in stars in a massive, central elliptical 
galaxy. For both the galaxy and the cluster systems the effect 
of cooling is only weakly dependent on radius and is significant 
even at r„; r . Note that axis ratios are not constant as a function of 
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FIG. 3. — Top and middle panels: Minor-to-major axis ratio, c/a, as a func- 
tion of radius for a cluster-size (top panel) and the galaxy-size (middle panel) 
halo. The thick solid lines correspond to the DM, while the thin solid lines 
show the combined c/a of DM and stars in simulations with cooling. The 
dashed lines show the profile of the DM in the adiabatic simulations. Bottom 
panel: The c/a profiles of merger remnants. We show remnants from four 
mergers: inclined halo+disk (solid line), coplanar halo+disk (short-dashed 
line), NFW halo (dotted line), and contracted halo (long-dashed line). Thick 
lines show c/a for DM only, and thin lines show c/a for stars in the cases with 
initial stellar disks in the inclined (solid line) and coplanar (short-dashed line) 
mergers. We show the profiles at a time ~ 8 crossing times of the remnant. 



radius and that changes are not mono tonic. Different regions of 
an object may be flattened to different degrees. This indicates 
that observational estimates of halo shapes in the disk region 
(e.g., by disk flaring, Oiling & Merrifield 2000) may be some- 
what different from results at larger distances. Note that despite 
the significant flattening of the baryons in the central 0.1r v i r of 
the galaxy simulation, the DM distribution around the disk is 
almost spherical. 

It is interesting to ask if the ellipsoid of the DM halo is 
aligned with that of the baryons. For clusters, in which most of 
the baryons in the center are in stars, the major axes of the stel- 
lar and DM distribution are approximately aligned at all radii. 
For example, the major axis of the central cluster galaxy is well 
aligned with the inner DM halo. However, we find that in clus- 
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ters that contain massive substructures in their outer regions, 
the direction of the major axis often changes dramatically at 
r > 0.5r vir . In the galaxy simulation, the minor axis of the DM 
distribution in the vicinity of the disk is aligned with the minor 
axis of the disk. However, the flattening of the DM distribution 
is small (Figure 3). Interestingly, at r > 0.2r v ; r , the direction of 
the major axis of the DM halo changes and is nearly perpendic- 
ular to the disk. 

We examined the evolution of the merger remnants and found 
that their shapes evolve in their outer regions for ~ 8 cross- 
ing times or ~ 14-18 Gyr. This indicates that the shapes of 
the outer regions of the cosmological halos are evolving at all 
epochs. The bottom panel of Figure 3 shows the axis ratios of 
the remnants in controlled merger experiments after 8 crossing 
times, when the evolution has ceased. The shapes of the DM 
halo merger remnants vary, depending strongly on the presence 
of a disk component and the relative inclination of disks prior 
to the merger. Mergers of halos with different central density 
profiles produce remnants with very similar axis ratios. When 
a disk component is present, the shape of the remnant DM 
halo depends sensitively on the initial relative inclination of the 
disks. Inclined disk mergers lead to a very spherical stellar rem- 
nant and a correspondingly more spherical DM halo, compared 
to the halo-only cases. Coplanar disk mergers lead to a very 
disklike stellar component (small c/a) and a DM halo that is 
nearly as flattened as in the halo-only mergers. We discuss the 
implications of these results in the next section. 

4. DISCUSSION AND CONCLUSIONS 

We show that halos in cosmological simulations with cooling 
are considerably more spherical than those in dissipationless 
simulations. The difference decreases with increasing radius 
but can be significant even at the virial radius. This is some- 
what surprising because cooling affects the mass distribution 
appreciably only in the inner <~ 10% of the virial radius. The 
condensation of baryons due to cooling leads to a more con- 
centrated distribution of DM, as it responds to the increasing 
gravitational field of baryons in the center (Blumenthal et al. 
1986). Thus, dissipation results in a significantly more cen- 
trally concentrated mass distribution and a deeper gravitational 
potential. Dubinski (1994) showed that this leads to the evolu- 
tion of the halo toward a more spherical shape in a few cross- 
ing times, arguing that as the central condensation grows, the 
overall potential becomes rounder. This shifts the boundaries 
between orbital families markedly, decreasing the fraction of 
regular box orbits that serve as the backbone of a triaxial mass 
distribution (Gerhard & Binney 1985; Udry & Martinet 1994; 
Barnes & Hernquist 1996; Merritt & Quinlan 1998; Valluri & 
Merritt 1998). 

In hierarchical models of structure formation, halos grow 
via a sequence of violent mergers and periods of slow accre- 
tion. Although cooling can gradually make a halo more spheri- 
cal, subsequent mergers can produce highly elongated remnants 
(e.g., Moore et al. 2003), erasing the effect of dissipation dis- 
cussed above. If no significant cooling occurs after the last ma- 
jor merger, 5 the triaxiality of the halo will be largely determined 
by the merger. Hence it is important to consider how cooling 
affects the shapes of merger remnants. To this end, we analyze 
a suite of controlled merger simulations of pure DM halos and 
halos with embedded disks. 

5 For example, if the merger occurs after most of the gas is converted to stars 
or the cooling time in the merger remnant is long. 



Cooling can directly affect the shape of stellar remnants dur- 
ing mergers (Barnes & Hernquist 1996). However, a large amount 
of cooling gas may be needed for this to significantly affect the 
shapes of DM halos. Indeed, we compare remnant shapes in 
mergers of disk+halo systems in which disks contain both stars 
and a modest amount of gas (10% of the total disk mass), with 
and without cooling. This comparison shows that the effect of 
dissipation on the shapes of dark halos during the mergers of 
these stellar-dominated disks is negligible. 

On the other hand, our merger simulations demonstrate that 
the shape of the remnant DM halo is affected by the presence of 
cold stellar disks. Figure 3 shows that the shape of the DM dis- 
tribution in the remnant is always correlated with the shape of 
the stellar distribution. Inclined disk mergers produce a nearly 
spherical stellar remnant, while coplanar disk mergers yield a 
very oblate stellar remnant with c/a ~ 0.3. The correspond- 
ing DM distribution is also more spherical in the inclined disk 
merger. The shape of the stellar remnant is thus quite sensitive 
to the mutual orientation of the merging disks. It is likely to be 
also sensitive to their internal properties. Indeed, in an inclined 
halo+disk merger with the same orbital parameters as before, 
but with a factor of 3 thicker stellar disks, the stellar remnant is 
less spherical and the halo axis ratios are lower by ~ 0. 15 in the 
inner r < 0.1r v ; r . If a spherical remnant is formed, the potential 
in the central region, where stars dominate gravitationally, will 
also be spherical. This likely destroys the regular box orbits of 
DM particles and drives the DM distribution to a more spherical 
configuration. This effect is not confined to the vicinity of the 
stellar remnant. The DM distribution may be affected at larger 
radii because DM particles in CDM halos are, on average, on 
eccentric orbits with a median apocenter-to-pericenter ratio of 
rapo/V = 6/1 (Ghigna et al. 1998). 

Our results indicate that the effect of cooling during quies- 
cent evolution and its indirect effect during late-time disk merg- 
ers are both important in determining the shape of the final mass 
distribution. These effects lead to halos that are systematically 
more spherical than those measured in dissipationless CDM 
simulations. This has several important implications for com- 
parisons with, and interpretations of, observations. The mean 
axis ratios of MW-size halos may be expected to be (c/a) ~ 
0.7-0.8 in their inner regions compared to (c/a) ~ 0.6-0.7 
estimated from dissipationless A^-body simulations. This can 
explain the observed coherence of the tidal streams of the Sagit- 
tarius dwarf spheroidal galaxy without requiring the MW halo 
to be an outlier in the distribution of halo shapes. We may also 
expect any intrinsic correlations between the shapes of neigh- 
boring galaxies to be mitigated. The effect of triaxiality on the 
average lensing cross section of clusters may not be as strong as 
previously thought (Dalai et al. 2003; Oguri & Keeton 2004). 
More spherical inner halos may influence the dynamics of bars 
and stars in galactic disks (e.g., El-Zant & Shlosman 2002; 
Gadotti & de Souza 2003) and may limit the mass of nuclear 
supermassive black holes by inhibiting their fueling (Merritt & 
Quinlan 1998). 

The effect of cooling is likely to be mass dependent because 
the efficiency of cooling and star formation, manifested in the 
average mass-to-light ratio M/L, is a fairly strong non-monotonic 
function of mass (e.g., van den Bosch et al. 2003). Current 
cosmological simulations that include gasdynamics suffer from 
the "overcooling problem": the fraction of baryons in stars and 
cold gas is at least a factor of 2 higher than observed for the 
systems of the mass range that we consider (Lin et al. 2003). 
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The effect on the axis ratios is therefore likely overestimated 
and substantial improvements in cosmological simulations are 
required to predict the shape distribution of CDM halos accu- 
rately. Nevertheless, the effect that we find is so strong that 
even with significantly reduced fractions of cold gas, we ex- 
pect the change in the axis ratios in the inner regions of halos 
to be > 0.1, comparable to or larger than the dispersion in the 
distribution of measured axis ratios in dissipationless N-body 
simulations. 
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